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Abstract 

We generalize the inverse patchy colloid model that was originally developed for heterogeneously 
charged particles with two identical polar patches and an oppositely charged equator to a model 
that can have a considerably richer surface pattern. Based on a Debye-Hiickel framework, we 
propose a coarse-grained description of the effective pair interactions that is applicable to particles 
with an arbitrary patch decoration. We demonstrate the versatility of this approach by applying 
it to models with (i) two differently charged and/or sized patches, and (ii) three, possibly different 
patches. 
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I. INTRODUCTION 


Colloids with patterned surfaces (commonly referred to as patchy particles) are considered 
as very versatile building entities whose shape and surface decoration can be designed such 
that they support the self-assembly of target structures with desired properties [H [2]. The 
self-organizing behavior of these units is based on their ability to form directional and 
highly selective bonds via the specific interaction features of the different surface areas. 
The versatility of this class of particles can be even enhanced if their surface decoration 
is characterized by charged regions, since in this case different surface areas can either 
be mutually attractive or mutually repulsive according to their respective charges. Many 
naturally occurring units at the nanoscale, such as proteins, virus capsids, and spotted 
vesicles, are characterized by differently charged moieties [3HZ3; in simple model systems 
an heterogeneously charged surface can be realized in the so-called inverse patchy colloids 
(IPCs) 0. 

IPGs are charged spherical colloids decorated by oppositely charged regions, realized, for 
instance, by polyelectrolyte stars adsorbed onto the particle surface pun]. Based on the 
features of the complex units presented in Ref. [9], a related model has been developed [H], 
where two identical, charged patches are located on the poles of an oppositely charged 
colloid. Closed expressions for the single particle potential and for the pair potential of 
two interacting IPCs can be written using a standard formalism of electrostatics within 
the Debye-Hiickel approach: the resulting analytical description of the pair interactions is 
based on approximating - under high screening conditions - the ratio between spherical 
Bessel functions of consecutive orders with a Yukawa-like term [8]. In order to make IPC 
systems amenable to theoretical approaches and/or computer simulations a coarse-graining 
description was also introduced, leading to considerably simplified expressions for the pair 
interaction between IPCs with two symmetric patches [8j. These interactions, that depend 
on both the relative distance and the relative orientation of the two IPCs, are able to form 
bonds in a highly selective manner, leading to complex self-assembly scenarios [HHE]. 

The rapid progress in synthesizing heterogeneously patterned colloids with desired patch 
decorations PZH23] urges to push forward the described IPC model towards a more sophis¬ 
ticated framework where particles can have an arbitrary patch decoration. The present 
contribution is dedicated to present a general formalism that allows to calculate the effec- 
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tive potential between two arbitrarily decorated IPCs. Similar to the original derivation for 
two identical patches [5], we replace the differently charged regions of the particle by point 
charges that are located in their respective centers of charge. For the ensuing charge dis¬ 
tribution we solve Maxwell’s equations, along with the boundary conditions at the surface, 
and obtain the potential created by a single IPC: expanding both the charge density and 
the ansatz for the potential in terms of spherical harmonics we can derive for a particular 
patch decoration a set of equations for the expansion coefficients. Based on this potential, 
we can then calculate the effective interaction between two IPCs. The algorithm, although 
amenable to an arbitrary patch decoration, becomes rapidly involved with an increasing 
number of patches and/or the complexity of the patch decoration. In this contribution we 
demonstrate the versatility of our approach by considering in detail the two-patch case with 
a size and/or charge disparity in the patches; moreover we extend the present approach to 
IPCs with a symmetric patch decoration of three identical patches. 

The manuscript is organized as follows. In section |TT| we provide the theoretical descrip¬ 
tion of the effective interactions between two IPCs with polar patches that are characterized 
by two possibly different sets of parameters. Following the Debye-Hiickel approach, we first 
- analytically - calculate the electrostatic potential around a single two-patch IPC (see 


subsection IIA) and then we - numerically - determine the effective potential between two 


IPCs (see subsection IIB). In section III we introduce the coarse-grained model and we 
motivate the choice of its parameters with the mapping procedure that links the analyti¬ 


cal to the coarse-grained description. In section IV, we consider a broad selection of IPC 


systems and we show results for the Debye-Hiickel potential (see subsection IV A) and for 


the coarse-grained potential (see subsection IVB). In subsection IVA, we distinguish three 
cases: IPCs with identical (i.e. symmetric) patches as in Refs. [TT, 12], IPCs with patches 
carrying different charges, and IPCs with patches having different surface extensions; re¬ 
sults for the case of symmetric patches are compared to those obtained with the Yukawa 


approximation used in Ref. [8] • In subsection IV B we consider the case of patches that are 
asymmetric in charge. Our conclusions are outlined in section [V] Finally, in the Appendix 
we demonstrate how to extend our approach to IPCs with three identical patches. 
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II. THE DEBYE-HUCKEL DESCRIPTION 


We consider an IPC of radius a surrounded by a dielectric solvent. According to Gauss’ 
law the electrostatic field outside a sphere with a homogeneously distributed surface charge 
is identical to the held generated by a point charge positioned in the center of the dielectric 
colloidal particle. When a colloid has a heterogeneously charged surface, the charge of 
the different areas (either the patches or the bare colloid) can be replaced by a discrete 
distribution of point charges positioned at the respective centers of charge inside the particle: 
while the center of charge of the colloid coincides with the center of the colloidal sphere, the 
centers of charge of the patches are located inside the particle in a well-defined geometry (that 
exactly reflects the patch decoration), each at a specific distance from the particle center. 
The inside of the colloid is inaccessible to both the co- and counter-ions of the electrolyte 
solution, while outside the particle these ions are part of the surrounding medium. For sake 
of simplicity we assume that the dielectric permittivity, e, has the same value both inside 
and outside the colloid. 


In Subsection |II A| we derive the screened electrostatic potential generated by a single 
particle with a discrete distribution of charges: the potential is calculated both inside and 
outside the colloid by expanding it in terms of spherical harmonics; then, by imposing 
electrostatic boundary conditions, a set of linear equations for the expansion coefficients is 
obtained. The system of equations can be numerically solved, leading eventually to the single 


particle potential. Based on this information, we numerically calculate in Subsection IIB 
the pair potential between two identical particles. 

The present approach can be applied to a colloidal particle decorated by an arbitrary 
number of patches, n p , at well-defined positions. As the algorithm becomes rapidly involved 
with increasing n p , we focus here on IPCs with two differently charged or sized patches. To 
demonstrate the versatility of the concept, we briefly outline the case of IPCs with three 
patches in the Appendix. 
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A. The electrostatic potential around a single two-patch IPC 

1. Calculation of the potential inside the colloid 

We describe the heterogeneously charged colloid as a dielectric sphere with a discrete 
charge distribution p(r, 9 , ip) located in its interior. For sake of simplicity we set the elemen¬ 
tary charge q e to unity. The potential generated by the distribution of the charges inside 
the particle satisfies Poisson’s equation (Gaussian units are used in the following) 

4tt 

A$ (1) (r,0,<^) = — p{r,Q,ip)\ (1) 

here the charge density, p(r,9,ip), is the result of the central charge of the colloid, Z c , and 
the two out-of-center charges of the patches, Z P1 and Z p2 , respectively; the latter ones are 
assumed to be positioned at distances cq and a 2 from the particle center, in directions 
opposite to each other. In spherical coordinates the charge density is thus 


p(r, 0, ip) = -Z C S( r) - Z V1 \s{r - 0^5(9 - ^)%>) - Z p2 \s(r - a 2 )S(9 - ^)S(<p - vr). (2) 

CL 2 Zd 6^2 Zd 

The general solution of Equation (JTJ) is the sum of the general solution of the corresponding 
homogeneous equation (also known as Laplace’s equation) and a particular solution of the 
inhomogeneous equation. The first contribution, $^ m (r, 0, </?), is 

OO +£ 
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the Ye m (9,ip) being the spherical harmonics. Since the potential inside the colloidal sphere 
should be well-behaved at r = 0, the coefficients Ci m are set to zero. 

A particular solution of Equation (jlj), < hp!^. t (r, 9,p), depends on the specific charge dis¬ 
tribution p(r, 9, ip). In our case - Equation (J2|) - we end up with the following expression 
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where the star denotes complex conjugation. We note that for two identical patches, the 
expansion of the potential in terms of spherical harmonics can be replaced by an expansion 
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in terms of Legendre polynomials [ 8 ]. Here the coordinate system is fixed such that one 
patch is located at ip = 0 and 9 = |. This system is used in all following steps. 

The solution inside the colloid, 9, ip), is finally given by 

$ ( 1 ) (r, 9 , p) = ^ m (r, 9, p) + d>part ( r , 9 , p). (5) 


2. Calculation of the potential outside the colloid 


Outside the colloid the dielectric solvent has to be taken into account: it consists of a 
large number of co- and counter-ions which will be treated within a mean held approach 
such that Boltzmann statistics apply. The characteristic features of the dielectric solvent 
enter only via the inverse Debye screening length k m- At low ion densities, a linearization 
of the Poisson-Boltzmann equation leads to the Helmholtz equation for the potential outside 
the colloid, <fh 2 )(r, 9, p), 

A<f>( 2 )(r, 9, ip) = n 2 & 2 \r, 9, p). ( 6 ) 


In its most common form this equation has a different sign and is solved by an expansion 
in terms of the spherical Bessel functions ji(x) and ye(x). Using the substitution k —> in in 
the expression for the standard solution, we arrive at the following solution for Equation ([ 6 ]) 

00 +i 
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where i is the complex unit and 
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with n = 0,1, 2, 3,.... Since $^ 2 ^(r, 9, <p) has to be regular for r —» 00 , while j n {i^ r ) and 
y n (inr) both diverge at large distances, we introduce - following standard procedures - 
suitable linear combinations of the two sets of functions, namely 
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which are commonly known as Hankel function of the first and second kind, respectively. 
Since the Hankel functions of the second kind do not decay exponentially with the distance, 
only the Hankel functions of the first kind are retained in the following. 

Thus the potential outside the colloid can be written as 

OO +£ 

9, <p) = EE Ci m hf\iKr)Ye m (d,tp). (12) 

£=0 m=—i 

3. Linking the two potentials via boundary conditions 


To solve the system of differential equations, we need to consider the proper boundary 
conditions. First the potential must vanish at infinity and has to be continuous at r = a. 
Second, in the absence of surface charges, the tangential component of the electrostatic field 
has to be continuous. Finally, the normal component of the displacement field must be 
continuous. Since the relative permittivity e is taken to be the same inside and outside of 
the colloid, these boundary conditions can be written as 


$ (1) (r, 9,ip)\ r=a = <f> (2) (r, 9,ip)\ r=a 
<9 r <F (1) (r, 9, ip)\ r=<T = d r & 2 \r,0,ip)\ r=a . (13) 

We now insert the expansions of $^^(r, 6 , ip) - see Equations (JdJ) and ([ 3 ]) - and of db 2 ) ( r — 
a, 9 , ip) - see Equation (12) - into these relations. Since the spherical harmonics are linearly 
independent, we obtain for each index combination (£, m) an equation for the corresponding, 
yet unknown coefficient that has to be solved numerically; the series expansions are truncated 
at a suitable upper limit, £ max . Once these coefficients are known, the resulting electrostatic 
potential around a single 1PC is thus available. Being interested only in the potential 
outside the colloid, we will omit in the following the superscript 2 and refer to the screened 
electrostatic potential in the region of interest as <f>(r, 9, ip). 

Since in our particular case the patches are arranged on the (6 = |)-plane, the 0-symmetry 
is conserved and only spherical harmonics with an even sum (£+m) contribute to the expan¬ 
sion. If, in addition, the patches are symmetric in size and charge, a further symmetry can 
be exploited; as a consequence, only functions with m mod n p = 0 need to be considered [8], 
where n p denotes the number of patches in the system. 

In Figure [I] we show the single particle potential <3>(r, 9, ip) generated by an IPC carrying 
two identical patches as a function of 9 , keeping r fixed to either r = a (left panel) or 
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r = cr + 4 /k (right panel). The system parameters used for the comparison are specified 
in the figure caption. In these panels we compare results calculated via the analytic route 
(presented here) and data obtained via the approximate approach put forward in Ref. [5]; 
the two different routes will be labeled henceforward as “numerical approach” (NA) and 
“Yukawa approximation” (YA), respectively. 

We note that even though the potential reported in Ref. |8j can formally be obtained 
by considering only the lowest-order spherical harmonics, these data cannot be reproduced 
within the present formalism by simply setting i max = 0. This discrepancy can be traced 
back to how the cut-off is conducted in the different approaches: in Ref. [Sj the sum over 
the Legendre polynomials is approximated after the total potential has been calculated 
analytically; in contrast, in the numerical calculation, the approximation needs to be carried 
out directly in the basis functions, before the boundary conditions are imposed. 



FIG. 1: Single particle potential $(r,0,cp) for an IPC with two identical patches; the potential is 
represented as a function of the angle 9 at a fixed distance: r = a (left panel) and r = cr + f (right 
panel); for two identical patches, the single particle potential does not depend on ip due to the 
azimuthal symmetry. The following set of parameters is chosen: Z c = —180, Z pi = Z p2 = Z p = 90, 
a\ = 02 = a = 0.6cr, kg = 5 and e = 80 (dielectric permittivity of water at room temperature). 
In both panels, the continuous red line corresponds to the numerical approach (NA) developed 
here - with cut-off l max = 80 for the number of spherical harmonics -, while the dotted green line 
reproduces the data obtained via the Yukawa approximation (YA) used in Ref. [Sj- 








B. The effective potential between two IPCs 


Once the single particle potential is known, the effective interaction energy between two 
identical colloids can be calculated using the same procedure as put forward in Ref. [ 8 j: 
we first determine the potential energy due to the presence of an IPC in the screened 
electrostatic field generated by another IPC (see Equations (5) and ( 6 ) of Ref. [ 8 ]) and 
then we symmetrize it; finally, the total interaction energy for a given particle-particle 
configuration is obtained as the average value over the two contributions divided by the 
absolute value of the minimum of the attraction. 

The described procedure is analytically applied in Ref. [ 8 ] due to the Yukawa approxi¬ 
mation operated on the single particle potential; however, in the present contribution we 
calculate the pair potential via a numerical route, thus retaining higher-order spherical har¬ 


monics. 


The radial and angular dependence of the resulting pair potential is investigated in Sec. IV 
on changing the microscopic parameters of the model; in particular, we focus on the effect 
of varying either the surface charges, i.e. Z c , Z Pl , and Z p2 , or the patch positions, i.e. cq 
and a 2 , while keeping the screening conditions, imposed by na, constant. 


III. THE COARSE-GRAINED DESCRIPTION 

We consider the coarse-grained description introduced in Ref. [ 8 ] for IPCs characterized 
by two identical patches and we extend this model to the case of two different patches. The 
coarse-grained description is designed to reproduce the same symmetries as the underlying 
microscopic system and is characterized by three independent sets of parameters: the inter¬ 
action ranges of the different surface regions, their interaction strengths, and their surface 
extents. 

The model features a hard spherical particle of size a carrying two interaction sites placed 
in opposite directions at distance a\ and <22 from the particle center; both distances are 
smaller than a so that the two sites are always located inside the colloid. As a consequence, 
the corresponding site interaction spheres (with radii pi and p 2 ) extend only partially out¬ 
side the hard core particle, defining in this way the polar patches; the respective surface 
extensions of the patches are characterized by the half opening angles 71 and 72 . 
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Since the characteristic interaction distances in the microscopic system are determined by 
the electrostatic screening of the surrounding solvent, all entities of the colloid are assumed 
to have the same interaction range, 8, irrespective of the surface regions involved in the 
interaction. The interaction sphere of the bare colloid has thus radius <7 + 5/2. 

By construction, the following relations between the geometric parameters of the model 
hold 


2 a 


cosy* = 


Pi 

o 2 + a 2 


a 


Pt 


2cr a i 


(14) 

(15) 


where z = 1,2; these relations are the generalizations of Equations (10) and (11) in Ref. [ 8 ], 
where only the special case cq = cq and p\ = p 2 was studied. By virtue of these constraints, 
the model is characterized only by two sets of geometric parameters: the particle interaction 
range 8 and the patch sizes 71 and 72 . Indeed, once the independent parameters of the 
model, 012 and P 72 , are defined, the physical parameters, 8 and 772 , are also fixed. 

In contrast, the energy parameters of the model are related to the charges involved in the 
interactions. While the center of the particle carries a charge Z c , the two sites carry charges 
Z Pl and Z P2 , respectively. These charges are responsible for the ratio between the attractive 
and repulsive contributions to the pair energy associated to the different (patch/patch, 
patch/bare and bare/bare) interactions. 

The specific form of the pair potential is based on the postulate that each of the afore¬ 
mentioned contributions can be factorized into an energy strength and a geometrical weight 
factor, the latter one being given by the distance dependent overlap volume of the involved 
interaction spheres. More specifically, beyond the hard core repulsion, the pair potential 
between two IPCs at distance r is given by | 8 j (omitting the arguments for brevity) 


U = 


4na 3 S ij U ij W ij 

0 


if 2 a < r < 2 cr + 8 
if r > 2 a + 8, 


(16) 


here i ( j ) specifies either a site or the center of the first (second) IPC, Wij is the overlap 
volume of the corresponding interaction spheres, and is the energy strength of the ij 
interaction. We note that, while the are constants fixed by the mapping procedure, the 
- as well as the potential U - depend on both the inter-particle distance and the relative 
orientation of the two IPCs. Although complex functions of the particle-particle distance 
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and orientation, the Wij can be written as simple functions of the site-site, site-center and 
center-center distances, as reported in Ref. [S]. 

Both the length and the energy scales involved in the coarse-grained model are imposed 
by the mapping schemes outlined in Ref. [8]. For what concerns the length scales, we assume 
that the coarse-grained interaction range is proportional to the Debye screening length ac¬ 
cording to the following relation: kS = n, where k is determined by the screening conditions 
hid = m; thus 5 = —a, where m and n are not necessarily integer numbers. For a given <5, 
the angular patch extents 7 i i2 are defined by the choice of either ai >2 or p 12 . The overlap 
volumes w^j for a given configuration are fixed by the choice of the geometric parameters. 
On the other hand, the energy scales are imposed by the mapping scheme referred to as 
“max” in Ref. [8]: the energy strengths Uij are fixed by considering characteristic reference 
configurations of two IPCs at contact and then imposing that the contact values of the 
Debye-Hiickel and the coarse-grained potentials are equal. For a given set of geometric pa¬ 
rameters, the overall particle charge Z tot = Z c + Z pi + Z p2 determines the set of values 
for the chosen reference configurations. 

For IPCs with two symmetric patches the three configurations depicted on the right hand 
side of Figure [2j namely the polar-polar (PP), the equatorial-polar (EP) and the equatorial- 
equatorial (EE) configuration, are enough to fully describe the inter-particle potential. In 
this case, the pair interaction energy is normalized with the value corresponding to the EP 
configuration mm- For the full description of IPCs with asymmetric patches (either in 
size or in charge), six characteristic particle configurations are needed: in particular, in the 
PP and EP configurations the existence of two different patches must be taken into account. 
In the following we denote differently charged or sized patches as patch 1 and 2, and we thus 
distinguish six characteristic configurations, i.e. PP(ll), PP(12), PP(22), EP(1), EP(2), 
and EE. In this case, the interaction energies are normalized by the most negative value 
occurring (i.e. the most attractive interaction), corresponding either to the EP(1) or the 
EP(2) configuration. 

IV. RESULTS 

In the following we evaluate the effective potentials between different pairs of two-patch 
IPCs in characteristic particle-particle configurations. More specifically, we investigate in 
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subsection IV A the radial and angular dependence of the Debye-Hiickel potential, while we 


discuss the coarse-grained counterpart in subsection IV B 


A. The Debye-Hiickel potential 


In the following we first consider the case of IPCs with symmetric patches and compare 
the Debye-Hiickel potentials obtained with the NA developed here to the results obtained 


within the YA proposed in Ref. [8j (see IV A 1); we then discuss the NA-based potentials 


for IPCs with asymmetric patches, either in charge (see subsection IVA2) or in size (see 


subsection IV A 3). 


1. Symmetric patches 

We first consider overall neutral IPCs with two identical patches, i.e., we set Z tot = 0, 
Z p = Z pi = Z p2 , and a = a 1 = a 2 . For this case we compare the radial and angular 
dependence of the pair potentials obtained within the NA developed in this contribution 
and within the YA put forward in Ref. [8], To this end we consider a selected set of 
parameters: in particular, we fix the screening conditions kct, the charges Z c and Z p , and 
the patch-center distance a such that we recover the same set of model parameters referred 
to as 45n in Ref. HU; this label identifies the selected model by the patch size of the related 
coarse-grained description, i.e., 7i = 72 = 45°, and by the overall particle charge, where “n” 
stands for neutral. In the following, we continue to use this labeling scheme. As shown in 
Fig. 0 the results obtained within the two approaches are sufficiently close to each other: 
only in the case of the PP repulsion we observe that it is slightly stronger within the NA; 
for the chosen set of parameters (specified in the caption of Fig. [ 2 ]), the NA contact energy 
of the PP configuration is by a factor « 1.15 bigger than the corresponding YA-value. 

As the contact energy, V r 2cr , for the different particle-particle configurations characterizes 
the pair potential, we calculate such a value on varying the charge imbalance between 
the patches and the bare colloid, starting from the 45n model and increasing the colloid 
charge Z c ; the charge imbalance is thus quantified by the overall particle charge Z tot . The 
corresponding data, labeled with 45, are reported in Fig. [3j in the main panel we show the 
NA-data, while in the inset YA-results are displayed. Since for the EP configuration VA is 
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FIG. 2: Radial (left) and angular (right) dependence of the interaction energy V (in units of k-^T) 
of two overall neutral IPGs with two identical patches. The specific choice of the parameter set 
guarantees that data obtained within the YA correspond to the model referred to as 45n in Ref. HU- 
Z c = —180, Z pi = Z p2 = Zp = 90, oi = <22 = a = 0.44a, k<t = 5, and e = 80 (dielectric permittivity 
of water at room temperature). Results obtained within the NA are represented by continuous 
lines, while YA-based data are depicted by dotted lines. In the left panel, the potentials are shown 
as functions of the inter-particle distance r between two IPCs in three different mutual orientations 
(as labeled): the polar-polar (PP), the equatorial-polar (EP), and the equatorial-equatorial (EE) 
configurations, depicted schematically in the right panel; in the right panel, the potentials are 
represented for two IPCs at contact, i.e. at a distance r = 2a, as functions of the rotation angle u 
around an axis perpendicular to the plane of the figure: after a rotation of 90° one particle-particle 
configuration transforms into one of the other particle arrangements. 

always equal to unity by construction, only the contact values for the repulsive contributions 
are reported. Both approaches show the same trend on increasing Z tot : while in the overall 
neutral system the PP repulsion is much stronger than the EE repulsion, in overall charged 
systems the PP repulsion decreases whereas the EE one increases up to the point where the 
EE interaction becomes more repulsive than the PP one. The comparison between the NA 
and the YA reveals that even though this trend is present in both sets of data, the values 
of Ztot at which the PP and the EE repulsions are equal differ between the two approaches: 
for the chosen set of parameters (specified in the caption of Fig. [3]), Ztat/Z p ~ -1.11 within 
the NA, while this value is approximately -0.44 within the YA, indicating that the described 
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trend of the characteristic repulsions is more pronounced in the latter potential description. 



FIG. 3: Contact energies (in units of k&T) of two IPCs at contact in the EE and PP configu¬ 
rations (as labeled) as functions of the charge imbalance Z tot /Z p ; the contact energy between two 
IPCs in the EP configuration is by construction equal to -1 and is thus not displayed. The main 
panel reports the NA-results, while the inset shows the corresponding data obtained within the YA. 
For fixed screening conditions and particle charges, three patch-center distances (corresponding to 
three different patch sizes in the related coarse-grained models) are chosen such that the overall 
neutral systems studied in Refs. mm are recovered. For overall neutral IPCs, the following set 
of parameters is chosen: Z c = —180, Z P1 = Z p2 = Z p = 90, a± = 0,2 = a = 0.64 a (labeled 30 to 
specify the coarse-grained patch size 7 = 30°), 0.44cr (labeled 45 to highlight the coarse-grained 
patch size 7 = 45°), and 0.32(7 (labeled 60 to specify the coarse-grained patch size 7 = 60°), na = 5 
and e = 80 (dielectric permittivity of water at room temperature); for overall charged IPCs, the 
patch charges are fixed to Z p = 90, while the colloid charge is changed from Z c = —180 by Zt 0 t- 

The radial and angular dependence of the pair potential between charged IPCs charac¬ 
terized by almost identical EE and PP contact energies is reported in Fig. [4j for the chosen 
set of parameters (specified in the caption of Fig. [I| we fix in our NA and YA calculations 
= —1.11 and Z to t/Zp = —0.44, respectively; the latter system corresponds to the 
set of model parameters that is referred to as 45c in Ref. EH; again, such a label iden¬ 
tifies the selected model via the patch size of the related coarse-grained description, i.e., 
71 = 72 = 45°, and by the overall particle charge, where “c” stays for charged. As in the 
neutral case, the trends in the radial and angular dependencies of the potentials are the 
same for both approaches, the NA providing data that are only slightly more repulsive than 
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the YA-based results (for both the EE and the PP configurations). 



FIG. 4: Radial (left) and angular (right) dependence of the interaction energy V (in units of 
ksT) between two overall charged IPCs with two identical patches. The specific choice of the 
parameter set guarantees that data obtained within the YA correspond to the model referred to 
as 45c in Ref. [llj: Z c = —280 (for the NA) or —220 (for the YA), Z Pl = Z P2 = Z p = 90, 
ai = a -2 = a = 0.44er, no = 5, and e = 80 (dielectric permittivity of water at room temperature). 
Results for the NA are represented by continuous lines, while data obtained with the YA approach 
are depicted by dotted lines. In the left panel, the potentials are reported as functions of the inter¬ 
particle distance r between two IPCs in the PP, EP, and EE configurations (as labeled); in the right 
panel, the potentials are represented for two IPCs at contact, i.e. at a distance r = 2a, as functions 
of the rotation angle uj around an axis perpendicular to the plane of the figure: after a rotation of 
90° one particle-particle configuration transforms into one of the other particle arrangements. 

Finally, we studied the behavior of the contact energy as a function of the charge imbal¬ 
ance for other IPC models, specified by a different choice of the parameter a. We note that 
at fixed interaction range, the choice of a uniquely determines the patch size 7 of the related 
coarse-grained model. In Fig. [3] we report VY-values for IPCs that are labeled by 30 (with 
a patch-center distance bigger than the one of the 45 model and such that 7 = 30°) and 60 
(with a patch-center distance smaller than the one of the 45 model and such that 7 = 60°), 
assuming EE and PP configurations. The choice of the parameters (specified in the caption 
of Fig. [3]) is such that for Z to t = 0 the YA reproduces the interactions of the overall neutral 
models studied in Ref. [IT] (referred to as 30n) and in Ref. [12] (referred to as 60n); accord¬ 
ingly, their overall charged counterparts investigated in Refs. mm are reproduced by the 
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YA potentials for Z to t/Z p = —1.11 (model referred to as 30c) and —0.22 (model referred 
to as 60c). Despite a scaling factor, the NA and YA provide very similar trends for the 
repulsive contact energies in the PP and the EE configurations as functions of the charge 
imbalance. Thus, we can speculate that also for models 30 and 60 the full radial and angular 
dependencies of the corresponding potentials will not show significant differences between 
the two approaches, as shown for model 45. We can conclude that both the NA and the YA 
provide consistent results for the description of the inter-particle interaction between IPCs. 


2. Asymmetric patches in charge 


We now consider IPCs with two polar patches of equal size but different charge. In 
particular, we focus on the 45n model investigated in the previous subsection and we change 
the charges of the two patches, while keeping the central charge of the colloid fixed. In order 


pi. 


is 


to maintain a vanishing overall charge of the particle, the charge of one patch (say Z } 
increased from Z p to Z p (l + q), while the charge of the other patch (say Z P2 ) is decreased 
from Z p to Z p (l — q). In Fig. [5] we report the resulting contact energy values as functions of 


q for the six characteristic particle-particle configurations specified in Sec. Ill As expected 
the three PP energies, as well as the two EP energies, are identical for q = 0 and they 
deviate from each other as soon as q increases. We note that since the V^-value for IPCs 
in the EP(1) configuration is used to normalize the pair potential, the corresponding curve 
in the (V 20 -, q) -plane is of course a horizontal line; on the other hand, the EP(2) value for 
V 2 a increases from -1 to slightly positive values as the charge of the second patch decreases, 
i.e. on increasing q. We also note that the contact value of the EE repulsive interaction 
is almost constant - it changes from slightly positive to slightly negative values - as the 
patch-free zones are unaffected by the change of the patch charges, at least as long as the 
IPC is overall neutral. Finally, as Z P2 increases and Z Pl decreases, the PP(ll) and PP(22) 
contact values rapidly increase and decrease, respectively; in contrast, the PP(12) contact 
value, being given by the overlap of one decreasing and one increasing patch charge, features 
a much slower decrease. In Fig. [5]there is a narrow g-range where the curves of the PP(12), 
PP(22), EE and EP(2) contact energies intersect (i.e., at q* ~ 0.66 for the particular set of 
parameters chosen). Beyond this point the PP(22) energy becomes attractive although the 
patches still have the same sign. This value depends on the system parameters: at fixed 
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interaction range, the q -value at which all the interactions - except the reference EP(1) one 
- have very similar contact energies increases as the patch size decreases. Nonetheless, since 
at this point the charge asymmetry is already quite large, the potential description must be 
considered unreliable. 



FIG. 5: Contact energies V^o (hi units of ksT) between two overall neutral IPCs in six different 
mutual orientations (as labeled) as functions of the charge asymmetry parameter q (defined in the 
text). When q = 0, the model corresponds to the 45n model whose parameters are given in the 
caption of Fig. [2] The models with differently charged patches (q > 0) are obtained in the following 
way: the central charge is kept constant to Z c = —180 and the two patch charges, initially set 
to the same value Z pi = Z P2 = Z p = 90, are either increased or decreased by a factor q, i.e. 
Z Pi = Z P (1 + q) and Z P2 = Z p (l — q)\ Z p sets the charge unit. The vertical black line highlights 
the value of q = 0.22 for which the radial and angular dependence of the pair potentials is shown 
in Figure [6} 

In Fig. [6] we show the full radial and angular dependencies of the inter-particle potential 
for a IPC model with asymmetrically charged polar patches, specified by q — 0.22. An 
interesting effect of the asymmetry is that the positions of the various minima and maxima 
in the angular dependence of the potentials are slightly shifted to larger or smaller angles 
with respect to the symmetric case. A detailed view of the potentials for u> ~ 90° in the 
right panel of Fig. [6] shows that, for example, the minimum of the EP(1) attraction is not 
located at u> = 90° but rather at 83° for the chosen set of parameters. The optimal particle- 
particle configuration is thus such that the less charged patch of one IPC is closer to the 
more charged patch of the other IPC: for the specific set of parameters, the former is shifted 
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by 7° closer to the latter with respect to the perfect EP configuration; as a consequence the 
distance between the centers of the more charged patches on the two different IPCs is bigger 
than their distance in the perfect EP configuration. We observe that, since different shifting 
angles are found for different fashions of rotational freedom, the optimal relative orientation 
may not be covered by the rotations studied here. This asymmetry effect is more pronounced 
as the charge asymmetry parameter q increases, so that only a small g-range can be studied 
with a reasonable accuracy within our framework. 

Investigations on models 30n and 60n (not shown here) have highlighted that for a given 
g-value, the shift of the minima and maxima in oj is more pronounced as the patch size 
increases. 



FIG. 6: Radial (left) and angular (right) dependence of the interaction energy V (in units of 
k-QT) between two overall neutral IPCs at contact carrying two differently charged patches. The 
asymmetry between the two patch charges is specified by q = 0.22 (highlighted by the black vertical 
line in Fig. [5j. In the left panel, the potentials are reported as functions of the inter-particle distance 
r between two IPCs in six characteristic particle-particle configurations (as labeled); in the right 
panel, the potentials are represented for two IPCs at contact, i.e. at a distance r = 2a, as functions 
of the rotation angle uj around an axis perpendicular to the plane of the figure: after a rotation of 
180° one particle-particle configuration transforms into another. In the inset of the right panel a 
magnified view of the potentials around u = 90° is shown. 
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3. Asymmetric patches in size 


Finally, we consider IPCs with two polar patches of equal charge but of different size. 
We start again from the 45n model and now change both the interaction range of one polar 
site, e.g. the site corresponding to patch 2, and its distance from the particle center in such 
a way that the quantity a 2 + pi remains constant. This procedure guarantees that in the 
corresponding coarse-grained system, the interaction range of patch 2 is the same as the 
interaction range of the whole particle, i.e. 6 is fixed. For a specific interaction range, there 
is thus only one choice for CZ2 for a given p 2l and vice versa; as a consequence there is only 
one the patch size y 2 associated to the choice of p 2 and a 2 or, equivalently, of p 2 and 5. 

In Fig. [7] we report the behavior of the contact energy values V 2a of the selected IPCs 
as functions of p 2 . The vertical line in the figure marks the size of patch 1 as well as the 
size of patch 2 in the 45n model: when p\ = p 2 the three different PP repulsions coincide, 
while both EP attractions are equal to -1. When p 2 decreases (i.e. when q 2 decreases), both 
repulsions PP(22) and PP(12) increase, the first one faster than the latter one: when p 2 is 
decreased by about 20%, the PP(22) repulsion has increased its value by a factor of three, 
while the PP(12) repulsion has only doubled. In both cases the increase of the repulsion 
is due to the increase of charge per area related to the decrease of patch size for a fixed 
patch charge. We note that also that the corresponding PP(ll) value is slightly affected by 
a change in the size of patch 2, but its increase amounts to less than 10% of its initial value; 
thus it can be considered to be essentially constant. The same consideration applies for 
the EE repulsion, which changes with the increase of the size of patch 2 by less than 20%. 
Since the surface charge increases on decreasing the patch size, also the EP(2) attraction 
increases: it has more than doubled when p 2 is decreased by about 20%; of course, the EP(1) 
attraction is constant by construction. The opposite trends are observed when p 2 increases 
(i.e. when y 2 increases); however, the changes in the contact energy values are in this case 
less pronounced. 

B. The coarse-grained potential 

We finally compare the potentials obtained via the Debye-Hiickel approach and via the 
coarse-graining procedure for the 45n system with a charge asymmetry parameter q = 0.22; 
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FIG. 7: Contact energies V^o- (in units of k-^T) between two overall neutral IPGs in six different 
mutual orientations (as labeled) as a function of the radius of one patch, p 2 , while the radius of the 
other patch is fixed to p\ = 0.76<r. The charges of the patches are equal. The black vertical line 
marks the set of parameters that corresponds to the 45n model; left to the line the size of patch 2 
decreases, while right to the line its extent increases. 


data are reported in Fig. [8] and show a good quantitative agreement between the Debye- 
Hiickel description and the coarse-grained approach. The radial dependence of the potentials 
is displayed in the left panel of the figure, while their angular dependence is shown in the right 
panel of the same figure. By definition, at r = 2a and oj = 0° (i.e. when the two IPCS form 
the reference configurations) the coarse-grained data coincide with the results obtained in the 
analytical calculations since the contact energy values are fixed by the mapping procedure. 
For larger distances (left panel) or angles (right panel) the coarse-grained data follow the 
analytical ones quite well; however, it must be noted that the angular dependence of the 
coarse-grained potential around u> = 90° shows small differences with respect to the Debye- 
Hiickcl data. Around the minima and maxima obtained via the analytical approach, the 
coarse-grained potential forms plateaus that are due to the saturation of the overlap volume 


between the interaction spheres. The effects described in subsection IV A 2 are thus not 
reproduced by the coarse-grained model: while in the analytical calculations there is an 
optimal bonding angle, in the coarse-grained description there is a range of angles around 
the optimal one. As long as this range is reasonably narrow, Monte Carlo simulations based 
on the coarse-grained potential will give accurate results at finite temperatures. Indeed, 
when studying ensembles of IPCs in simulations under these conditions small deviations 
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from the analytical pair energy will affect the system behavior only to a minor extent, due 
to the role played by entropy. 



FIG. 8: Theoretical and coarse-grained pair potentials V (in units of k^T) for overall neutral IPGs 
with patches of size 45 and charge asymmetry q = 0.22. Results for the Debye-Hiickel potential 
are represented by dotted lines, while those for the coarse-grained interaction (for n5 = 2) are 
depicted by continuous lines. In the left panel, the potentials are reported as functions of the 
inter-particle distance r between two IPCs in six different mutual orientations (as labeled); in the 
right panel, the potentials are represented for two IPCs at contact, i.e. at a distance r = 2cr, as 
functions of the rotation angle oj around an axis perpendicular to the plane of the figure: after 
a rotation of 180° one particle-particle configuration transforms into another. In the inset of the 
right panel a magnified view of the potentials around cj = 90 is shown. 


V. CONCLUSION 

In the present contribution, we have proposed a nearly analytic method to calculate 
the effective interaction potential between inverse patchy colloids (IPCs) with rich surface 
patterns. Our description follows the Debye-Hiickel approach to determine the potential 
generated by one single colloid dispersed in a liquid solvent: the potentials inside and outside 
the colloid are calculated separately as expansions in terms of spherical harmonics and 
are then linked together by electrostatic boundary conditions; the resulting set of linear 
equations for the yet undetermined expansion coefficients is solved numerically. Once the 
single particle potential is known, the effective interaction energy between two IPCs is 
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determined as in Ref. [8]. The main advantage of not solving the equations analytically 
is that they can be kept as general as possible so that they can be easily applied to a 
wide range of different IPCs. The Debye-Hiickel potential is subsequently used to provide 
a suitably designed coarse-grained model with parameters that are directly related to the 
physical quantities of the underlying microscopic system, such as the screening conditions, 
the surface charges and the extent of the different surface regions on the colloids. The 
coarse-grained description is advantageous in many body simulations because the model pair 
energy between two particles is the sum over simple products of geometric (orientation- and 
distance-dependent) factors and energy factors, associated to the different ways of interaction 
that characterize the IPC. 

In this contribution, the pair potentials for neutral as well as overall charged colloids with 
two symmetric or asymmetric patches in size or/and charge are calculated. A comparison 
between the characteristic energies of the Debye-Hiickel description and of the coarse grained 
model shows good agreement. In an effort to demonstrate the versatility of our approach, 
the framework is also applied to the case of three possibly different patches. 

The proposed method will allow future studies on the self-assembly of heterogeneously 
charged colloidal systems, that are expected to be observed in experiments. 
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VI. APPENDIX: IPCS WITH 3 PATCHES 

To demonstrate the versatility of our approach we consider a particle of charge Z c and 
diameter 2a, which is decorated by three patches of charges Z Pl , Z P2 , and Z P3 . The patches 
are located on the equatorial plane of the colloid, their respective centers-of-charges being 
separated by the distances cq, a 2 , and a 3 from the center of the colloid; connecting lines 
between a patch center and the particle center enclose throughout 120°. 
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The source term of Equation (|TJ) is given by 


p(r,9,<p) = ~Z C 6( r) - Z Pl \s(r - ai)S(6 - ^)S(<p) 

CL -1 Z 


( 17 ) 


1 


7I\ 


- Z^Sir - ai)6(e -‘-)6 y) 


Zp, l ;,Yr - njl'M" - pi (p - 4 ( T j 


leading to the following particular solution of the differential equation inside the colloid 
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Combining expression (18) with the homogeneous solution given by expression (|3|, we obtain 
the potential inside the colloid < mb(r, 0,<f>). 

Outside the colloids, the solution of the Helmholtz equation (Jb]), <f>( 2 )(r, 9, 0), is given by 


expression (12) 


Imposing on 6, (/)) and $part(r, 6,<f>) the boundary conditions (|13|) at r = a leads 


,(2) 


to a set of linear equations for the unknown expansion coefficients that is solved numerically. 
In this way, we obtain the potential created by one three-patch IPC, $(r, 9, 0). It is worth 
noticing that for this patch decoration the axial rational symmetry is lost. 

In Figure [9] we display <f>(r, d, 0) for a IPCs with three identical patches (left panel) or 
three differently charged patches (right panel) at r = cr. 

The pair potential of two interacting IPCs is obtained by following the same procedure 
applied in the two patch case; the resulting interaction potential depends in a complex way 
on the center-to-center distance of the colloids and on their relative respective orientations. 

In parallel to the derivation of the Debye-Hiickel potential, we design the coarse-grained 
model along the same lines traced in the body of the paper for the two patch case: we 


make an ansatz for U, similar to the one in Equation (16), where the sum is now extended 


over all the possible pairings of the four interaction spheres. The evaluation of the overlap 
volumes Wij is straightforward given the relative position and orientation of the two IPCs, 
while the yet undetermined energy parameters u VJ are again fixed by considering well-defined 
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FIG. 9: The single particle potential, $(?’, 9,cj)), for an IPC with three identical (left panel) or 
differently charged (right panel) patches is reported at distance r = a from the particle center as 
a function of 9 and <j>. For identical patches, the chosen parameters are: Z c = —180, Z Pl = Z P2 = 
Z P3 = 90, a± = ci 2 = a3 = 0.44fj, ko = 5, and e = 80 (dielectric permittivity of water at room 
temperature). In the asymmetric case, the charges of the patches are chosen to be Z pi = 110, 
Z P2 = 70, and Z P3 = 90, while the other parameters are unchanged. 


reference configurations. When the three patches are identical (in size and charge), then 
three 'Ujj-values must be determined, corresponding to the pure bare/bare, bare/patch, and 
patch/patch interactions. The most characteristic configurations needed to evaluate these 


pure interactions are depicted in Figure 10 (labeled again EE, EP, and PP, respectively). 
As the choice of the reference configurations is to some extent arbitrary, it is equally valid 
to consider rotations of one of the two IPCs that preserve the pure overlap between the 
selected interaction spheres. The rotation axis for the chosen configurations is depicted in 
Figure [TT)} The average value of the analytic contact energies of all equivalent configurations 
is then used as the reference for the mapping to the coarse-grained contact energies. 

It must be noted that the mapping procedure can be performed also when the three 
patches are different. In this case, the w^-values to be evaluated adds up to ten and thus ten 
reference configurations with the pure EE, EP and PP interactions must be considered: one 
EE configuration, three EP configurations and six PP configurations. Due to the reduced 
symmetry as compared to the identical patch case, averaging over all possible rotations of 


one IPC with respect to the other is advisable. The same rotation axis shown in Figure 10 
can be used. 
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Side View 


Top View 



FIG. 10: Reference configurations for the mapping between the analytic and the coarse-grained 
pair potential of IPCs with three identical patches: characteristic configurations for the pure EE, 
EP and PP interactions are reported from top to bottom, as labeled; side views on the left, top 
views on the right. The second IPC can rotate with respect to the first one - preserving pure EE, 
EP and PP interactions - around an axis depicted as a black arrow in the side view or as a black 
dot in the top view. 
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